Quantum States and Phases in Driven Open Quantum Systems with Cold Atoms 



00 

o 
o 

(N 



S. Diehli'2, A. MicheRi^^ A. Kantiani-^, B. Kraus^^^ H.P. Biichler^ and P. ZoUeri-^ 
^ Institute for Theoretical Physics, University of Innsbruck, Technikerstr. 25, A-6020 Innsbruck, Austria 
^Institute for Quantum Optics and Quantum Information of the Austrian Academy of Sciences, A-6020 Innsbruck, Austria and 
'^Institute for Theoretical Physics III, University of Stuttgart, Pfaffenwaldring 57, 70550 Stuttgart, Germany 

(Dated: March 10, 2008) 

An open quantum system, whose time evolution is governed by a master equation, can be driven 
into a given pure quantum state by an appropriate design of the system-reservoir coupling. This 
points out a route towards preparing many body states and non-equilibrium quantum phases by 
quantum reservoir engineering. Here we discuss in detail the example of a driven dissipative Bose 
Emstem Condensate of bosons and of paired fermions, where atoms in an optical lattice are coupled 
to a bath of Bogohubov excitations via the atomic current representing local dissipation. In the 
absence of interactions the lattice gas is driven into a pure state with long range order. Weak 
interactions lead to a weakly mixed state, which in 3D can be understood as a depletion of the 
condensate, and in ID and 2D exhibits properties reminiscent of a Luttinger liquid or a Kosterlitz- 
Thouless critical phase at finite temperature, with the role of the "finite temperature" played by 
the interactions. 
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I. INTRODUCTION 

In condensed matter physics quantum phases and as- 
sociated strongly correlated many body states are typ- 
ically prepared by cooling the system to low tempera- 
tures, where its properties are dominated by the ground 
state of a Hamiltonian, H \G) = Eg\G), i.e. by con- 
sidering a thermodynamic equilibrium situation, where 
p ^ _^ fQj. temperature T ^ 0. 

In particular, in the context of ultracold atomic quan- 
tum gases, much of the interest of the last few years 
has focused on engineering specific Hamiltonians based 
on control of microscopic system parameters via exter- 
nal fields [U, [3, H, 0, H, 0], opening the door to a 
quantum simulation of strongly correlated ground states 

0,ii0,[iil,[I3,El. 

In contrast, quantum optics typically considers driven 
open quantum system, where a system of interest is 
driven by an external field and coupled to an environ- 
ment inducing a non- equilibrium dynamics, with time 
evolution described by a master equation for the reduced 
system density operator (see e.g. |14|). 

(1) 
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Here H is the Hamiltonian of the driven system, while the 
Liouvillian C in Lindblad form represents the dissipative 
terms. The quantum jump operators C£ are system opera- 
tors as they appear in the interaction Hamiltonian for the 
coupling to the bath of harmonic oscillators, and describe 
the time evolution (quantum jump) of the system asso- 
ciated with the emission of a quantum into the harmonic 
oscillator bath with rate Kg in channel £. The validity 
of the master equation is based on the Born-Markov ap- 
proximation with system-bath coupling in rotating wave 
approximation, which in quantum optics is an excellent 
approximation because the (optical) system frequencies 
are much larger than the decay rates. For long times, the 



system described by ([T]) will approach a dynamical steady 
state, p{t) pss, which in general will be a mixed state. 
However, under special circumstances the steady state 
can be a pure state, pss = \D) {D\, where in the language 
of quantum optics \D) is called a dark state. Sufficient 
conditions for the existence of a unique dark state are 
(i) V£ ci\D) = 0, i.e. the dark state is an eigenstate of 
the set of quantum jump operators with zero eigenvalue, 
which will be compatible with the dynamics induced by 
the Hamiltonian if (ii) H\D) = E\D). Uniqueness of the 
dark state is guaranteed if there is no other subspace of 
the system Hilbert space which is invariant under the ac- 
tion of the operators ci [l^[l3|. In fact, it can be shown 
that for any given pure state there will be a master equa- 
tion so that this state becomes the unique steady state 

Here, we are interested in this novel possibility of 
quantum state engineering by designing jump operators 
such that one drives the system into a desired many- 
body quantum state. This non- equilibrium approach is in 
strong contrast to conventional Hamiltonian engineering 
methods, as standard thermodynamics concepts are not 
valid in this driven system, and the dynamics goverened 
by the master equation ([T]) is the only remaining prin- 
ciple determining the final state. While in quantum op- 
tics we know several examples of preparing single particle 
pure states dissipation, including dark state laser cool- 
ing to subrecoil temperatures [17, 18], it is of interest to 
extend these ideas to many body systems, dissipatively 
driving the system into entangled states of interest, or 
preparing non-equilibrium quantum phases in condensed 
matter systems. Furthermore, for the example of a dissi- 
pative driven Bose Einstein condensate (BEC) discussed 
below, but also for stabilizer states in a system of spins- 
1/2 or qubits living on a lattice (iGj . the dissipation can 
be chosen to be quasi-local, i.e. the jump operators act 
non-trivially only on a small neighborhood of particles. 

As an illustration of a many particle dark state, we 
discuss below a dissipatively driven BEC, where for non- 
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interacting atoms a pure state exhibiting long range order 
is generated as the steady state by quasi-local coupling to 
an environment with finite correlation length. Applying 
standard linearization schemes in the weakly interact- 
ing situations, allows us to determine the solution of the 
master equation ([T]), and we find that the steady state 
exhibits similar properties as bosons in thermal contact 
to a heat bath: in 3D the effect of the interaction can 
be understood in terms of a depletion of the condensate, 
while in f D and 2D the system exhibits properties remi- 
niscent of a Luttinger liquid or a Kostcrlitz-Thouless crit- 
ical phase at finite temperature. In particular, we give 
a physical realization in terms of cold bosonic atoms in 
an optical lattice by immersion in a superfluid bath [l^ . 
As a second example we discuss a master equation whose 
steady state corresponds to an ?7-condensate [20,,^], i.e. 
a state of long range order of paired fermions, which re- 
markably corresponds to an exact excited eigenstate of a 
Hubbard Hamiltonian with repulsive or attractive inter- 
actions in d dimensions. 



II. LONG RANGE ORDER BY LOCAL 
DISSIPATION 

A. Dissipative Driven Condensate 

Let us consider the dynamics of N bosonic atoms on d- 
dimensional lattice with M'^ lattice sites, lattice vectors 
ex and spacing a (see Fig. la). We assume that the 
coherent motion of the atoms can be described a single 
band Hubbard model with Hamiltonian 

H = Ho + V^-jJ2alaj + ^uY,'^r'^l (2) 

where Hq represents the kinetic energy of bosons hop- 
ping between adjacent lattice sites with amplitude J, 
and V is the onsite interactions with strength U, and 
cti {a]) are bosonic destruction (creation) operators for 
atoms at site i. A physical realization of this situation 
is achieved by loading cold bosonic atoms into an optical 
lattice. By cooling to temperatures T ^ the system 
is prepared in the groundstate, which for noninteracting 
atoms is |BEC) = flq^Q |vac) /\/iVT which corresponds to 
a state with macroscopic occupation of the quasimomen- 
tum q = 0. Here Oq — J2j aje^^^^ / Vm^ is the destruc- 
tion operator for quasimomentum q in the Bloch band. 
For weak interactions we have in 3D a BEC, and a quasi- 
condensate in ID and 2D, while with increasing interac- 
tions the system can undergo a quantum phase transition 
to a Mott phase [11,0, ll. 

In contrast, we are interested here in a dissipative Hub- 
bard dynamics modelled by a master equation ([1]). In 
atomic physics a bath will typically couple to the atoms 
via the atomic density, Ui = ajat, as in the case of deco- 
herence due to spontaneous emission in an optical lattice, 
or for coUisional interactions. This will tend to dephase 



the condensate, and can heat the system. In contrast, our 
goal is to couple the system to a bath so that the system 
is driven to a pure many body state by quasi-local dissi- 
pation. This is achieved, for example, by chosing jump 
operators 

ce = Cij = (^al + a]^ (a^ - a^) (3) 

acting between a pair of adjacent lattice sites £ = {i,j) 
with a dissipative rate = k It is then easy to see that 
the state |BEC) satisfies (i) V(i, j) (a, - a/) |BEC) = 0, 
and (ii) is an eigenstate of the kinetic energy Hq |BEC) = 
Ne^^o |BEC), where Cq = 2J^^ sin^ qeA/2 is the single 
particle Bloch energy for quasimomentum q. In other 
words, I BEC) is a many body dark state corresponding 
to a state of long range order for non-interacting bosons 
on the lattice. In fact (see Appendix A), this state is 
the unique steady state of this master equation, and any 
initial mixed state will evolve for long times into |BEC). 

The key in obtaining a state of long range order as the 
steady state is to couple to the bath involving the atomic 
current operator between two adjacent lattice sites; the 
concept that dissipative coupling to the current operator 
stabilizes superconductivity is well knwon in condensed 
matter [2^ . The jump operator Cij describes a pump- 
ing process where the second factor — aj annihilates 
the anti-symmetric (out-phase) superposition on the pair 
of sites {i,j), while aJ + a] recycles the atoms into the 
symmetric (in-phase) state. Loosely speaking, we can in- 
terpret this process as a dissipative locking of the atomic 
phases of two adjacent lattice sites, which in turn results 
in a global phase locking, i.e. a condensate. Note that 
in view of (a^ — aj) |BEC) = it is the destruction part 
of the jump operator which makes |BEC) the dark state, 
while any linear combination aJ and of recycling op- 
erators will do, except for a hermitian which would 
lead to a pure dephasing but no pumping into the dark 
state. 

The above discussion becomes particularly clear in a 
momentum representation suggested by the translation 
invariance. The dissipative part of the master equation 
takes on the form ([T]) with jump operators 

Q = Cq,, = ^(1 + e'('^-^)-)(l - e-^^^^)al_^a^ 

and Kg = K, which makes the appearance of the dark state 
q = decoupled from dissipation particulary apparent. 

B. Implementation 

The above master equation can be realized by immers- 
ing atoms a moving in an optical lattice in a large BEC 
of atoms b [l^ . The condensate interacts in the form of a 
contact potential with interspecies scattering length aat 
with the atoms a, and acts as a bath of Bogoliubov excita- 
tions. This coupling provides an efficient mechanism for 
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FIG. 1: Driven Dissipative Condensate: (a) A lattice gas a 
is immersed in a condensate b, which acts on the links 
(crosses) of neighboring lattice sites i and j in the form a dis- 
sipative current, (b) Schematic realization of the effective dis- 
sipativ process in an optical super-lattice, which provides for 
excited states gapped by e and localized on the links of neigh- 
boring lattice sites (ij): A Raman laser couples the ground- 
and excited bands with effective Rabi-frequency fl and de- 
tuning A = u — e from the inter-band transition. The spatial 
modulation of the Raman-laser yields in a position dependent 
coupling, which excites only the anti-symetric component of 
atoms sitting on neighboring lattice sites i and j into the up- 
per band. The inter-band decay with a rate F back to their 
lower band is obtain via the emission of phonons into the 
surrounding BEC. 



decay of atoms a from an excited to lower Bloch bands by 
emission of a Bogoliubov quasiparticle. Atoms a moving 
in the lowest Bloch band can be driven by laser induced 
Raman processes to the first excited Bloch band, and can 
decay back to the ground band. This situation is rem- 
iniscent of optical pumping in quantum optics, or laser 
cooling [3, [13, [3 ■ There a laser excites electronic states 
of an atom, which return to the ground state by sponta- 
neous emission of a photon. This formal analogy suggests 
and justifies a description of a driven dissipative Hubbard 
dynamics in terms of a master equation. Thus our goal 
is to identify an optical lattice configuration which after 
adiabatic elimination of the excited Bloch bands in the 
limit of weak driving results in a master equation of the 
type discussed in Sec. Ill Al 

We consider a lattice as illustrated in Fig. la, with an 
additional auxiliary lattice site on each of the links. The 
optical lattice corresponding to a single link is shown in 
Fig. lb. It has the form of a A-system with the two 
Wannier functions of lattice sites 1 and 2 representing 
two ground states, and the auxiliary state in the middle 
representing an excited state. We drive this three-level 
system by Raman transitions from the two ground to the 
excited states with Rabi frequencies fl and —fl, respec- 
tively. The excited atom can decay back to 1 or 2 by 
emission of a Bogolibov quasiparticle. As is well-known 
from quantum optics, such a A-configuration supports a 
dark state: in the example of a single atom only the an- 
tisymmetric (out-of phase) state {a[ — a\) |vac) is excited 
by the laser, so that the atom is eventually "pumped" 
into the dark symmetric (in-phase) state {a\ + a\) |vac). 
In general, laser excitation followed by return of an atom 



to site a = 1,2 will involve operators ajj(ai — 02), and 
- as shown in Appendix B - results in a Liouvillian with 
the structure 



a/3 



alia, - aj), p(4 " al')a/3 + h.c.(4) 



The coefficients C are related to the correlation function 
of the Bogoliubov reservoir, and in particular exhibit the 
correlation of emitted quasiparticles at lattice sites a and 
(3 as reflected by the correlation length of the reservoir. 
For wavelength Ab of Bogoliubov excitations larger or 
smaller than the optical lattice spacing a, spontaneous 
emission is either correlated or uncorrelated. However 
complicated C, the existence of a dark state is guaranteed 
by (ai — ttj) I BEC) = 0, a property which follows from 
the laser excitation step. We will confine our discussion 
below to the master equation with jump operators ([3]) 
with qualitatively similar results expected for the general 
case. 



III. COMPETITION OF HAMILTONIAN AND 
LIOUVILLIAN DYNAMICS 

For a realistic system with a finite interaction V, see 
Eq. ([2]), the |BEC) is no longer a dark state of the master 
equation: the interactions tend to localize the particles, 
while the dissipative terms tend to enforce a pure con- 
densate. In general, the competition between these two 
incompatible dynamics results in a mixed state. Below 
we present linearized theories, which allow us to solve for 
the density matrix p, as well as to study the correlations 
functions under the time evolution. 



A. Mean field theory 

For weak interactions, one can expect that the pure 
I BEC) state is only weakly perturbed with the zero mo- 
mentum mode Oq^o still macroscopically occupied. We 
follow therefore the standard Bogoliubov prescription 
and replace in the master equation the zero momentum 
mode by its mean value, i.e., = \/ uqM'^ with np the 
condensate density {n — uq n). In leading order, the 
jump operators in the master equation ([T]) reduce to 



a]. 



(5) 



and the coupling rates ki = 16nK sin^ (qeA/2). Here, 
the index i — {q, ct} with cr = ±1 characterizes the bosonic 
operators at = Oq cr = (aq+(Ta_q) /\/2. In addition, 
the Hamiltonian in Eq. ^ simplifies to 



H — <{ei + U n) a\ai 



Un 



{4 



(6) 



with 

tion over the index £ avoids these double countings. It 



Eq. Note that aq 



(Ta_q^CT and the summa- 
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follows that the master equation decouples for each mode 
ai and is quadratic in the bosonic operators. The anal- 
ogous master equation in quantum optics is well known 
from parametric amplification. The general solution is 
given by a mixed Gaussian state, which in steady state 
takes the form p = Z^^Ilepe with 



pi = exp 



(7) 



with describing a finite occupation of the squeezed 
operators bi = e""^* cosh(6'f) ai + e"*" smhiOi) a\. The 
squeezing parameter 9^ and (3^ are given by the relation 

cosh^ {29,) = coth^ (/3,/2) = + (8) 



with El — y/eJ+2Une£ = Eq the Bogoliubov energy. 
The phase cot = (ei+nlJ) j Kg plays a minor role. From 
these relations, we recover the pure |BEC) in the limit 
of vanishing interactions Unj k ^ 0. On the other hand 
for fixed interaction and dissipation k, the modes hg es- 
sentially reduce to the well known Bogoliubov modes for 
small momentum Iqe^l ^ ^/U J / k, while the parameter 
takes the form (3g ~ Eg /T^a with an effective temperature 



= Un/2. 



(9) 



The density matrix for the low momentum modes is 
therefore indistinguishable from the thermal state of a 
weakly interacting Bosc gas with the role of a "finite tem- 
perature" played by the interactions Un. 

This similarity of the driven system with a thermal Bo- 
goliubov state naturally rises the question on the validity 
of the mean field approximation: the ansatz assumes a 
small depletion tt-d of the zero momentum mode. Using 
the above solution for the steady state, the condensate 
depletion takes the from 



no 



dq {Unf 



(10) 



with va the volume of the Brillouin zone. This expression 
strongly depends on the dimension of the system: in 3D 
the depletion remains finite and small for weak interac- 
tions Un/J <C 1. On the other hand, in one- and two- 
dimensions an infrared divergence appears indicating the 
absence of a macroscopic occupation and the breakdown 
of mean-field theory. Consequently, we find that a dis- 
sipatively driven system exhibits the same behavior as a 
bosonic system in themal contact with a heat bath, where 
the appearance of long-range order at finite temperature 
is only possible above the lower critical dimension d=2. 

The solution to the master equation also allows us to 
study the relaxation into the steady state, see Appendix 
C. The build up of the macroscopic occupation no(t) 
obeys the long time behavior 



The condensate no(t) approaches the steady state accord- 
ing to a power law. In the strict absence of interactions 
the behavior is modified to np— rio(i) ^ i^'^l'^ . The slower 
approach to equilibrium for the interacting system results 
from the scrambling of the particles via the interaction. 



B. Lower dimensions d = 1,2 

In lower dimensions, the interaction drives strong 
phase fiuctuations, which destroy the macroscopic occu- 
pation of the condensate. However, these fluctuations 
are only relevant on distances larger than the coherence 
length ^ J jUna [a denotes the lattice spacing). Conse- 
quently, the formation of a local condensate still takes 
place on shorter distances, while only the phase between 
these local condensates is destroyed by the fluctations. 
The influence of these phase fluctuations can be stud- 
ied in a long wave length description by introducing a 
smoothly varyirig phase field and density field Ui with 
\(\)i,ni\ = [24| . This behavior of the dissipatively 
driven system is in close analogy to the thermodynamics 
of interacting bosons giving rise to Kosterlitz-Thouless 
critical phases in 2D, and Luttinger liquids in ID. The 
jump operators simplify to 

Cij = (n^ - rij) - 2m {(^^ - 4>j) , (12) 

while the Hamiltonian ([2]) reduces to the harmonic model 



(13) 



no - no{t) 



(11) 



Consequently, the master equation becomes again 
quadratic and introducing new bosonic operators dg with 
£ — {q, a} allows us to decouple the master equation 
for each mode £: we define the bosonic operator dg — 
{dq + ad-q) with dq = (riq/y^ — iV2n0q) /\/2 . 
Surprisingly, with this definition the jump operator re- 
duces to cg = dg and the master equation ([1]) in the long 
wave length limit |q|a < ^Unj J is mapped to the same 
form as the master equation in the previous section, with 
the operator ag replaced by the new operator dg and the 
chemical potential Vn replaced by Un — eg/2. Conse- 
quently, we find again the exact solution to the master 
equation, which allows us to characterize the state via 
its correlation functions, see Appendix C; such a charac- 
terization of a states in cold gases has attracted a lot of 
interest recently 25, 26, 27, 28]. 

First, we analyze the steady state. The interesting 
correlation function in lower dimensions is the Green's 
function 

G^,{x,t) = {a,{h)a]{to)) (exp[i {^,{ti) - <t>,{to))]) 

with X the distance between the lattice sites i and j and 
t — \ti — t{)\. In the long wave length limit with x, tc^ 
■y/ J /Una, the smooth part of the correlation function 
recovers rotational and translation invariance; here c = 
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FIG. 2: Appearance of quasi long-range order during the time 
evolution: the correlation function Gt{x,0) is shown for var- 
ious times ten/2 = 0, 10, 10^ 10^, 10*, 10^, 10*^, oo. The in- 
tial disordered state has a correlation length ^ = 2a, and 
the system parameters are chosen at Teff/(4rKT) ~ 1/18 and 
xo = 0.55a. 



\/2U nJa is the sound velocity. Then, simple analytic 
results are found in 2D (with the help of the quantum 
regression theorem) giving rise to quasi long-range order 

a,(x,o) ^ I ^»^* (14) 

I (ro/i)*^KT x-^ct 

vifith the Kostcrlitz-Thouless temperature T^t = ttJu ^ 
Teff = Un/2, and the short distances scales xq ^ c/kh 
and To ~ (Kna/c)'^a/c. On the other hand, in one- 
dimension we find an exponential decay both in space 
and time with coherence length = AcK/nT^a, with 
the parameter K = iTy^2Jn/U playing the role of the 
Luttingcr parameter. The concept of the effective tem- 
perature is particularly efficient in the low dimensional 
systems, where the correlations are dominated by the low 
energy phase modes, as implied by the exponential (ID) 
and algebraic decay (2D). Noticeably, the functional de- 
pendence of the correlation function is determined by the 
ratio U/ J alone, while the dissipative coupling strength n 
only gives rise to non-universal prefactors. This behav- 
ior is a consequence of the low momentum behavior of 
the dissipative damping Kq ^ versus the linear sound 
spectrum of the Hamiltonian (|13p . 

Finally, we study the time evolution of the spatial 
correlation function in two-dimensions i.e., Gt(x,0) = 
(oio])* with the average defined by the density matrix 
p{t). The system is initially prepared in a disordered 
state characterized by a small correlation length ^, giv- 
ing rise to an exponential decay of the correlations. The 
exact time evolution then allows us to study the appear- 
ance of the quasi-long range order Eq. (|14p from this ini- 
tially uncorrelated state. For t ^ x"^ / (IGku), we obtain 



in leading order 

Gt(a;,0) - (x/xo)"^ e-^'/^l (15) 

The length scale within which the system exhibits quasi- 
long range order is Xt = 2{Tr^'^ unat)^^^ and increases in 
time with the universal power xt ~ t^/*, see Fig. [51 



IV. CONDENSATE OF INTERACTING 
FERMIONIC DOUBLONS 

As a second example, we consider the dissipative 
preparation of an ry-state, which is an exact excited 
eigenstate of the d-dimensional two-species fermionic 
Hubbard-Hamiltonian [13] 

Hfh^~J E /U. + C^E4/u/'i/«T- (16) 

This state is created by the 7y-operator, 77^ — 
J2i 'f'iVi /M'^/'^ with the local doublon operators ri\ = 

fiT-^h' '^^'^ 'pi — denoting a sign alternating between 
sites in a checkerboard pattern. The A^-ry-state is created 
by A^-fold application of 77^ , yielding an excited eigenstate 
oiHpH, {vT |0),with energy NU. This excited state is 
a condensate of doublons into the quasimomentum-state 
(tt, . . . , tt) in the corner of the d-dimensional Brillouin- 
zone. It exhibits superfluidity, with non-decaying off- 
diagonal long-range order in any spatial dimension. To- 
gether with the operator J2a flafiiy^ properties of rj^ 
allow the construction of symmetry generators of HpH, 
which has been used to investigate high-Tc superconduc- 
tivity (c.f. [21]). 

Following the lines presented above, we construct 
quasi-local operators at each pair of lattice sites having 
the A^-77- state as a dark-state. As is easily verified, the 
quasilocal operators Ci = c\j^ {k — 1,2) given by 

4'^ = U-rl''j){v^+m)^ (17) 
= n,^flfn+n,^flM, (18) 

fulfill these requirements. Small-scale numerical sim- 
ulations for open-boundary systems then show that a 
Liouvillian ([1]) constructed from these quantum jump- 
operators and the Hamiltonian Hfh drive any initial 
p(0) into the A^-77-state, assuming that N-^ — at all 
times. The result may be interpreted in the quantum 
jump picture: HpH generates configurations with spin- 
up and spin-down particles on adjacent sites from any 
initial configuration. These configurations are then cap- 
tured by cj^- , associating them into doublons. Subse- 
quent action of c^j' then generates the desired 7/-state by 
phase-locking, just as in the case of the jump-operators 
for the BEG. 
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V. CONCLUSIONS AND OUTLOOK 

We have discussed a scenario where many body quan- 
tum states and entangled states are prepared as dark 
states in a non-equihbrium driven dissipative dynamics 
with quasi-lcoal dissipation. While the present work has 
focused on condensed matter aspects of realizing non- 
equlibirium (quasi-)condensates of interacting bosons 
and paired fermions in optical lattices, the present ideas 
are readily extended to spin systems, and promise a new 
avenue towards preparing interesting entangled states of 
qubits for quantum information [l6|. On the atomic 
physics side, control via external fields offers interesting 
new possibilities of engineering a broad class of quantum 
jump operators ([31), one example being phase imprinting 
tti Oie"^' with a laser. 
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APPENDIX A: CRITERION FOR THE 
UNIQUENESS OF THE STAT IONARY STATE 
(SEE SEC. ITTA)) [16] 

The BEC is a dark state for the set of jump operators 
ce- For the uniqueness of this solution we further need it 
to be the only dark state. For a fixed particle number N, 
the first factor in ci (see Sec. (Ill A[ ) has no eigenvalues 
and in particular, no zero eigenvalues. Thus, in order to 
identify dark states \D) with zero eigenvalue, wc may re- 
strict to the equation (ui — aj)\D) = V(i, j). Taking the 
Fourier transform, this translates to (1 — e*''''^°)aq|D) = 
Vq. Thus the BEC state with q = is the only dark 
state. 

For the dark states to be the only stationary solution of 
p — Cp we need to show now that if there exists another 
stationary state, then there must exist a subspace of the 
Hilbert space 7i, which is invariant under any application 
of the operators c^. Equivalently, in |16j] we show the 
following 

Theorem 1. Let D he the space of dark states, i.e. 
cgD = 0. // there exists no subspace 5 C 7i with S J- D 
such that ciS C S Vq, then the only stationary states are 
the dark .states. 

In order to show that the jump operators Eq. (??) 
indeed lead to the unique BEC steady state solu- 
tion we construct a polynomial operator 0{ci) with 
{BEC\ 0$) 7^ for all 1$) e S" = \BEC). For this 
piupose we work in the momentum space representation 
such that 

with e^ characterizing the different lattice vectors con- 
necting nearest neighbour sites. A basis in the Hilbert 



space can be written as |{nq}) = nq('^q)"'' There- 
fore the polynomial operator O — Y[q Cq'' provides a fi- 
nite overlap with the BEC. A general state can be written 
as |<i>) = X]/{nq} K'^q})- Choose one basis state with 
nonvanishing coefficient with highest occupation in the 
zero mode. Then the polynomial operator O for this ba- 
sis state provides also a finite overlap of |$) with the dark 
state. 



APPENDIX B: IMPLEMENTATION OF THE 
BEC LIOUVILLIAN 

Coherently driven two-band Hubbard model - We out- 
line the mechanism for engineering Liouvillians driving 
into a BEC. We focus on a one-dimensional system, 
where the system atoms a are moving in an optical su- 
perlattice Vi^pt — Y!,n=i^nS>^'0-'^{nTTz/a)/n'^ with lattice 
spacing a along the direction z, while being tightly con- 
fined by a harmonic potential with oscillator frequency 
iu± in the transverse directions, x and y. The lattice 
depths V2 > Vi > of the superlattice are c hosen such 
that the vibrational spacings Sw„=o,i ~ 2 •\/^T^^2±^4) 
about the individual wells at positions zqj = ja and 
zi j = (j + l/2)a are much larger than their gap e ~ 
Vi + h{uji - wo)/2, cf. Figmb). Here E, = h'^n^ /2maa'^ 
denotes the lattice recoil and rua the mass of the atoms 
a. The low-energy dynamics of the atoms a is then given 
by a two-band Bose- Hubbard model 

n,q {cti} 

where al^ ^ — X^jli ^"^■'"^l j/v'Af creates a Bloch-wave 
with quasi-momentum < q < ir/a and in the lower and 
upper band, n = 0,1, respectively, and aj^ ^ creates an 
atom a in the site at position z„ j. We denote the disper- 
sion relations for the two bands by en,q — —^Jn cos{qa) + 
ne with J„ the nearest neighbor inter-band hopping rates 
in band n and e the band-separation. The second term 
in Ha with at = {ni,ji) account for intra- and inter- 
band interaction of two atoms a in the superlattice in 
terms of the respective on- and off-site shift C/{ai}- The 
dominant intra- and intra-band interactions are given by 
density-density interactions ^ j Unal^ju"^ j/2 and ~ 

^0.1 '^0 j'A j'^'i-j'^Oj, respectively. Here, C/„ denotes 

state-dependent on-site shifts, and Vn the intra-band off- 
site shift for nearest neighbors (ij), cf. « — j = 0, 1. 

A key ingredient for the BEC Liouvillian is a selective 
coherent drive between the two bands of the Hubbard 
model, which couples the antisymmetric superposition 
on a pair of lower sites to the upper level on the link in 
between. This is achieved via a Raman laser setup in the 
form of a set of standing waves locked at the position of 
one of its minima zij. This allows to realize an effective 
dynamical coupling Vias(i) = hil cos{Tr z / a) sin(ujt) with 
the same periodicity as the optical lattice. Here is the 
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two-photon Rabi frequencies, and we denote the Raman- 
detuning from the upper band by A — lu — e/h. For a 
weak driving field, <ti £/fi, the Raman-drive results in 
effective inter-band coupling 

3 

=hne~^-' ^(1 ~ e^n<iU+^/a^o,, + H.c. 



where we made use of the rotating wave approxima- 
tion and dropped AC-Stark shifts, since they only renor- 
malize the gap ~ Vi of the original superlattice, and 
fj — ftf is the effective Rabi-frequency up to a Franck- 
Condon factor /, which in terms of the Wannier-functions 
Wn{z — Zn,j) for the band n localized at z = Znj reads 
f = J dzw{){z) cos{TTz/a)wi{z — a/2). The alternating 
signs in the interband-couping iJias for atoms sitting on 
neighboring sites in the lower band n — originates from 
the modulation of the Raman-coupling. This results in 
a excitation of asymmetric (and in particular antisymet- 
ric) superposition of atoms sitting on adjacent sites from 
the lower band to the upper band, while their symmet- 
ric superposition is dark with respect to Ha, leading to 
the BEC dark state. From momentum representation of 
Hias(t) we see that the two counter-propagating compo- 
nents e^''^^/° provide a inter-band coupling of modes dif- 
fering by 7r/a, and their interference results in the Bloch- 
waves to be selectively excited from the lower to the up- 
per band based on their quasi-momentum q. In particular 
for the lower-band n = the two components destruc- 
tively (constructively) interfere for q = 0, (q = vr/a), 
while the upper-band displays the opposite behavior. 

Coupling to the phonon reservoir - The dissipative step 
is implemented by immersing the system into a large ho- 
mogeneous 3D condensate of a distinguishable species of 
atoms b, which acts as a reservoir of phonon modes . 
The corresponding phonon Hamiltonian is given in the 
Bogoliubov approximation by Hi, — J^k^o -E'fc&k^k, where 
&k creates a Bogoliubov excitation with momentum k 
and energy Ek = [{hkcb)"^ + {h^k'^ /2mbYY^'^ with mt the 
mass of the atoms b and Cb their speed of sound in the 
condensate. The atoms a and b interact via a contact- 
interaction with a coupling constant 17^6 = 27r S^OabZ/iafc 
given in terms of their intra-species scattering length Oat 
and their reduced mass /Zab = mairi}, / (rna+mb) . Expand- 
ing these density- density interaction in terms of the fluc- 
tuations about the condensate wave-function we obtain 
to first order (apart from an overall state-independent 
mean-field shift gabNaPb proportional to the condensate 
density pb) an effective coupling of the atoms a to the 
Bogoliubov excitations in the form 

Hab = gab ^ V PbSkAlb]^ + H.c. 
k/O 

Here 5*;^ = fi^k'^ /2mbEk is the static structure factor of 
the BEC and Aj^ is the displacement-operator with mo- 
mentum k for the atoms a associated with the recoil from 



the emission of a Bogoliubov excitation into the conden- 
sate. 

In the following we are interested in the effective dy- 
namics of the system atoms a and consider the BEC as a 
reservoir of Bogoliubov excitation at essentially zero tem- 
perature, since under typical experimental conditions one 
can achieve temperatures Tb <C e/fcb- We integrate out 
the bath dynamics in the Born-Markov approximation, 
and obtain a Master-equation for the density operator 
p{t) of the atoms a (within the rotating wave and in- 
dependent rate approximation with respect to the laser 
excitation) , 

^ ^^[Ha{t),p] + 1^ |||2([Ak,p4(£;,)] + H.C.), 

k 



with Ha{t) = Ha + H\iis{t) and where ^^(e) de- 
notes the Fourier-component of A]^ with frequency e/h 
with the frequency Ek/h given by the dispersion rela- 
tion of the Bogoliubov excitations. Given that intra- 
band dissipative processes are suppressed by momen- 
tum conservation p^ . we focus on the inter-band de- 
cay, and within the rotating frame approximation write 
the Master-equation with = ^ G 



k,g '-'l^q+fe^'^Cg 

and their Fourier-components given by A\J^E) « 
Al[T:5{Elh-£/h) + ir/{E/h-£/h)], where we ne- 
glected corrections Ja,Ua ^ £ the spectra and in the 
following will drop the Lamb-shifts, as they amount to 
(small) second-order shifts. 

We take the continuum limit of k and exploiting the 
radial symmetry of the Bogoliubov spectrum E^, perform 
the integration over k are left with a integral over its 
azimutal angle, which we rewrite in terms of kz/k as 



cp=y:-^ 



H.C. 



assuming a tight transverse confinement, e 
hu!± = h^ /2maa\, we approximate the Bloch- function 
by a Gaussian of width a±_ in transverse direc- 
tion times periodic component along z, $o-,g(r) ~ 
'/'(T.g(-z)e~''^/^°^ /(tto^)^/^ and perform the polar integra- 
tion in G'-^'"' as 



G 



(1,0) 
k,g 

(1,0: 



ft'o- /dz0r,+,^(z)e^'=^^0o,,(^) 



Long wavelength (Super /Sub-radiant) limit - In gen- 
eral the recoil along z for the bloch-wavefunction has to 
be computed numerically. However in the deeply bound 
limit a rough estimate may obtained by taking (orthogo- 
nalized) gaussians of width for the wannier functions 
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and focus on nearest neighbor contributions, which yields 



f(l,0) 
Jk,q 



e-'"^" I dze'^''-wi{z)wo{z - {21 - l)a/2) 

z=o,i 

2exp[-(A:aj2)2] sin2(A:a/8)e-^'?"/2 



sinh[(a/4a^)2]/ cos(fca/4 - qa/2) 
Thus the master-equation in momentum space reads 

/ dfc.e-'=^('^^-'^i)/2sin4(^)x 



(Y_^)cos(Y-¥) 



X 



H.c, 



withKo = 27rfcogfco^"''""'j-^^/?i^wosinh^[(a/4a2)2]. We no- 
tice that for ko < tt/o, one obtains a collective intra-band 
decay decribed by a set of jump operators Ck with mo- 
mentum transfer < fco < n/a along the lattice as 

with Kfc = Koe-'='("'-"'i)/''sin2(fca/8). 

Short wavelength limit - For /cq > '^7r/2 the sum over k 
spans over n Brillouin zones, resulting in construct and 
destructive interference, which one associates to the fi- 
nite correlation length of the Bogoliubov modes in the 
bath. In particular, for fci 3> fco = n/a we notice that 
the summation in the master-equation runs over several 
Brillouin zone, which effectively suppresses the non-local 
/ long-range decay. Thus it is convenient to rewrite the 
master-equation in position space, which yields the decay 
of the general form 



H.c 



i,j £,£' 



where the summations and £,£' are over the lattice 
sites and links, respectively, and in the limit k —^ oo the 
decay rates take the simple form 



To / dzwo^jWi,ewi,e'Woj', 



where Tq = QabSki/^ ~ S'^^ in terms of the being the 
static structure factor Ski momentum ki and Wn.i = 
Wn{z — Zn^i) denote the Wannier function for the lower 
(n = 0) and upper (n = 1) band, respectively. We notice 
that they are translationally invariant and rapidly fall 
off with increasing separation of {jii'j')- Thus we can 
restrict ourself to the largest quasi-local contributions, an 
in particular the ones for neighboring lattice sites yield 



j p=0,l 



70 
71 



H.c. 



The first term describes the uncorrelated decay from the 
link j to the left site j and to the right site j + 1 with 
^ ^fj ~ ^j+ij+i' "^^^ second one accounts for the 
finite correlation between the two processes (cf. ~ Jo) 
and the third is an analogous one that accounts for cor- 
relations in the decay between different links to the same 
site, i.e. links j — 1 and j to site j (cf. ~ Ji). It is conve- 
nient to rewrite the three-terms in terms of symmetric / 
antisymmetric combinations of particle localized on ad- 
jacent sites in the two bands, cf. (a„jipi ±a„j+izpi)/\/2. 
This results in a master-equation in the form of a "chain" 
of concatenated A-systems and inverted A systems, 



^p = E 

3,± 



±7o) 



H.c. 



ai 



V2 



V2 



P<j — 



V2 

_i ± a 



V2 



corresponding to a master-equation with jump operators 
C^^'"'^ on links (A-systems) and vertices (V-systems) re- 
spectively, as 



V7 ± 271 



O'l.jiaij-i ±aij) 



In the independent rate approximation we thus obtain 
the master-equation for the set of laser-driven excitations 
and quasi-local decay (that is in the regime where the 
wavelength of the Bogoliubov excitations is much smaller 
than the lattice-spacing, fci ^ tt/o) by 



dp 1 



[Ha + iJias(t), p{t)] + Cap + CvP, 



dt ih 

C^p^ ^ [2Ci'^piCl")^-{iCl^n'Cl^\p}], 

where aij denote the summation over links and (con- 
nected) vertices for a = A and a = V, respectively. 

From the master-equation we see that (i) all excited 
states decay back to the ground-state and (ii) for U = 
the phase-locked ground state (BEC) is a dark-state, i.e. 
an eigenstate of the hamiltonian that is dark with respect 
to the excitation and the decay. 

Since all excited states are decaying, we adiabatically 
eliminate the excited band. In the limit that one has a 
weak far detuned laser, A ^ J7, [/„, J„, 7^, 7^, results to 
lowest order aij w ^^('^Oj ~ ao.j+i)/\/2A and thus the 
jump-opeators transform to (up to global phases) 



C 



'^iKj ± aoj+i)(aoj - aoj+i)/2. 



A(j.±) 

Cv(j,±) =\/ «±Ooj(-ao,i-i + (1 ± l)aoj T ai,j+i)/2. 



where we dropped the band-index n — and we intro- 
duced the rates — 'y'^{fl/A)'^, and remark that the 
latter correspond to the diagonalized form of Eq. ^ . 
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APPENDIX C: MOMENTUM SPACE 
CORRELATIONS 

In the linearized theory aU information is encoded in 
the first and second moments. Our construction imphes 
the vanishing of the first moments, (a|) = (a^) = 0. We 
may therefore concentrate on the time evolution equa- 
tions of the second moments, whose solution then al- 
lows to reconstruct the Gaussian density operator and 
gives access to the correlation functions. In principle the 
relevant single particle sector of the density matrix is 
mapped out by the correlations (aja^'), {aiair),h.c., but 
for our purposes it is sufficient to focus on the diago- 
nal entries. The corresponding equations arc obtained 
from the linearized version of Eq. ^ (cf. Sect. IIIip 



using the commutation relation [a£,a],] = ^o-.u' (^q-q' + 
crSq^^q'),[ae,ai>] = [a|,a|,] = = (q,o-)), with the 
result 

dt{alae) = 2iC/n((a^a^) - (a^aj)) - 4K£(aJa£), (CI) 
dt{aiai) = -4i[Un{{alae) + ^) + (ei + Un)(aiai)] 
-AK.t{aiag) 

where we omit the equation for {a\a\) trivially obtained 
form the second line. The equations of motion may be 
solved via the Laplace transform and form the basis for 
the computation of the density operator and the spatial 
and temporal correlation functions. 
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